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Abstract 

There is a gap between single-species model predictions, and empirical studies, 
regarding the effect of habitat fragmentation per se, i.e., a process involving the 
breaking apart of habitat without loss of habitat. Empirical works indicate that 
fragmentation can have positive as well as negative effects, whereas, traditionally, 
single-species models predict a negative effect of fragmentation. Within the class 
of reaction-diffusion models, studies almost unanimously predict such a detrimen- 
tal effect. In this paper, considering a single-species reaction-diffusion model with 
a removal - or similarly harvesting - term, in two dimensions, we find both posi- 
tive and negative effects of fragmentation of the reserves, i.e. the protected regions 
where no removal occurs. Fragmented reserves lead to higher population sizes for 
time-constant removal terms. On the other hand, when the removal term is propor- 
tional to the population density, higher population sizes are obtained on aggregated 
reserves, but maximum yields are attained on fragmented configurations, and for 
intermediate harvesting intensities. 

Key words: fragmentation, single-species model, reaction-diffusion, harvesting, 
spatial patterns, conservation biology 



1 Introduction 

The analysis of the effects of environmental fragmentation and variability on 
population densities and biodiversity has stimulated the development of many 
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spatially-explicit population models. In the modeling literatur e, posit i ve ef- 



fects of environmental variability have bee n recorded (se e e.g. iBolkerl . 12003 



Berestycki et al.l . |2005aj) . On the other hand iFahrigl (120031 ) , in a thorough bib- 
liography analysis, pointed out that most single-species modeling approaches 
lead to comparable conclusions regarding the detrimental effects of fragmenta- 
tion per se, i.e., a process involving the breaking apart of hab itat without loss 
of habitat. She noted that unlike the effects of habitat loss (see Saunders et all 



199ll . for a discussion of the consequences of fragmentation with habitat loss), 
and in contrast to current theory, empirical studies suggest that the effects 
of fragmentation per se are at least as likely positive as negative. The aim 
of this note is to make steps towards a reconciliation between the theory 
and empirical works on the effects of fragmentation, within the framework of 
reaction-diffusion models. 



React ion- diffusion models (hereafter RD models), although they sometimes 
bear on simplistic assumptions su ch as infinite v elocity assumption and com- 
pletely random motion of animals (jHolmesl . ll993l ). are not in disagreement with 
certain dispersal properties of populations observed i n natural as well as exper- 
imen t al ecolog i cal sys t ems, at least qualitatively (seelShigesada and Kawasakil . 



19971 : iTurchinl . Il998l : iMurravl . 12002 : lOkubo and Levinl . l2002f ). Furthermore 



these models often provide a good framework for rigorous investigation of the- 
oretical questions and derivation of qualitative as well as numerical results on 
population dynamics. In that respect, the effects of environmental fragmen- 
tation have been addressed in many theoretical studies based on such models 
over the last decades. 



Within the class of RD models with heterogenous coefficients, numerous re- 
cent works have emphasized the detrimental effect of environmental fragmen- 
tation per se on species persistence and spreading, in agreement with the 
other theoretical tendencies noted by Fahrig. In all these RD models, the pop- 
ulation growth rate function at a location x, r(x), was not constant, taking 
higher values in favorable regions than in unfavorable ones. Depending on the 
spatial arrangements of these regions, the modeled populations were shown 
to tend to extinction or survive, and to disperse at different speeds. In the 
particular case of 1-dimensional binar y environments (i.e. for r t aking val- 



ues in a set constituted o f two v alues) , ICantrell and Cosnerl (119891 . 120031 ) and 
Shigesada and Kawasakil (119971 ) have establish ed that concentrating a ll the 
habitat in a single patch improved persistence. iBerestycki et al.l (j2005al ) gen- 
eralized these analytical results to the iV-dimensional case , with more general 
growth rate functions. More recently iRoques and Stoical (120071 ) carried out 
more precise results regarding the negative correlation between persistence 
and fr agmentation. Habitat fragmentation has also been shown, firs t num eri- 
cally (IKinezaki et all 120031 ). and then analytically (1E1 Smaili et all 120091 ) to 
negatively affect population spreading speed. 
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Nevertheless, to our knowledge, there do not exist detailed demonstration of 
positive effects related to fragmentation derived f r om R D models in two space- 
dimensions, although in a related work, iNeubertl (120031 ) shows positive effects 
of reserve fragmentation from RD models with harvesting terms, in the sense 
that fragmented reserves sometimes maximize the yield. However, in his work, 
which is carried out in a one-dimensional space, no explicit conservation of the 
of the reserve is assumed. 



In this paper, we study single-species RD models in two-space dimensions with 
a spatially-homogeneous growth term, and spatially heterogeneous removal 
terms — Y which can be, for instance, interpreted as harvesting terms. The re- 
gions where Y = hence correspond to protected regions or similarly reserves. 
In our models, Y = Y (x, u) depends on the location x, and can depend on the 
population density u at this location. Typically, when Y does not depend on 
the population density, it corresponds to a constant-yield harvesting strategy. 
In this constant number of individuals are re moved per unit of time . 

This is the case when a quota is set on the harvesters (|Robinson and Redford , 



199ll ; iRobinson and Bodmerl . Il999i ; IStephens et all 120021 ) . Even in the absence 



of such imposed quotas, harvesters often increase their effort to maintain a con- 
stant yield. A good example is provided b y the high-trophic l evel fi shes catches 
data in the North Atlantic reviewed by IChristensen et al.l (120031 ). which de- 
scribe a decline of the biomass of one half from 1950 to 1990, while the catch 
remains the same. The function Y can also be taken to be proportional to the 
density u, corresponding to a proportional harvesting strategy. Then, locally, 
a constant proportion of the population is removed per unit of time, corre- 
sponding to a constant effort of the harve sters. These two harves t ing st rategies 
have been investigated by the authors (IRoques and Chekrounl . 120071 ). in in- 
homogeneous environments. It was shown, through analytical and numerical 
analysis, that aggregated habitat configurations gave better chances of popu- 
lation persistence, respecting the tendency found in the modeling literature; 
see also lOruganti et al.l (120021 ) for other mathematical results on these models. 



We conduce here a different analysis. Not only we investigate the effects of the 
spatial arrangement of the harvesting term rather than those of the growth 
function, but we also focus on other quantities than simple persistence. Consid- 
ering protected regions with a fixed total area but with gradually fragmented 
shapes, we analyze the intertwined effects of fragmentation and harvesting 
intensity on both the population size and the quantity of harvested indi- 
vi duals. To do tha t we u se the stochastic model of landscape generation of 
Rogues and Stoical (120071 ). that is in complete agreement with the concept of 
fragmentation per se, and we show that fragmentation of protected regions can 
in fact be beneficial to the modeled population and to the harvesters. These 
results demonstrate that two-dimensional RD models with harvesting terms 
can support a dual effect of fragmentation per se, positive as well as negative. 
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Typically, it is shown, for instance, that fragmented reserves lead to higher 
population sizes for time- const ant removal terms. On the other hand, when 
the removal term is proportional to the population density, higher population 
sizes are obtained on aggregated reserves, but maximum yields are attained 
on fragmented configurations, and for intermediate harvesting intensities. 



2 Materials and methods 



2.1 The model 



The idea of modeling population dynamics with reaction-diffusion models has 
begun to develop at the beginning of the 20th century, with rando m walk 
theorie s of o rganis ms, introduced by Pearson and Blakemanl (119061 ). Then, 



Fisherl (119371 ) and iKolmogorov et al.l (119371 ) used a react ion- diffusion equa- 



tion wit h homo geneous coefficients as a model for population genetics. Later, 
Skellaml (Il95ll ) examined this type of model, and he succeeded to propose 
quantitative explanations of observations for the spread of muskrats through- 
out Europe at the beginning of 20th Century. Since then, these models have 
been widely used to explain spatial propagation or spreading of biological 
species (bacter ia, epidemiological agents , insec t s, fishes, mammal, plant s , etc. 



see the books IShigesada and Kawasaki! . 119971 ; iTurchinl . 1 1998c iMurrayl . 12002 
Okubo and Levinl . 12002k for review) . 



Ignoring age or stage structures as well as delay mechanisms or Allee effects, 
the classical Fisher-Kolmogorov model, in two space-dimensions, can be writ- 
ten as follows: 



du 
~dt 



- DS7 2 u = ru(l - u/K), t > 0, xGftC 



where u = -u(t,x) corresponds to the population density at time t and po- 
sition x = (xi,x 2 ). The left-hand side of ([[]) corresponds to the diffusion 
equation, and simply describes the redistribution of organisms following un- 
correlated random walks where V 2 stands for the spatial dispersion operator 
d 2 u d 2 u 

+ — — — . The diffusion coefficient D measures the individuals rate 



V 2 u 



dxl 



-1 

of movement, r > is the intrinsic growth rate of the population and K > 
corresponds to the environment carrying capacity. 



The domain Q is considered bounded, and we assume reflecting boundary 
conditions: 



du . 



0, for x e dSl, 
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where dQ is the domain's boundary and n = n(x) corresponds to the outward 
normal to this boundary. Thus, some part of the boundary can be considered 
as an absolute barrier that the individuals do not cross, like coasts, and other 
parts of the boundary can be seen as regions where as much individuals exit 
the domain as individuals enter the domain. 

At this stage no environmental fragmentation is present in this class of models. 
To introduce it, we adopt a perturbative approach which consists in subtract- 
ing a spatially-dependent term to the right-hand side of equation ([I]): 

— - DV 2 u = ru(l - u/K) - F(x, u), t > 0, xGllC R 2 . (2) 



It could seem artificial to keep constant the parameters r, D and K whereas 
the removal term Y is spatially-dependent. The ecological interpretation of 
such a framework is no more than the consideration of heterogeneously dis- 
tributed harvesting in homogeneous media, a classical set-up in fisheries for 
instance. In fact, even in the case of a population living in a heterogeneous en- 
vironment, the objective of such assumptions is to separate the effects in order 
to facilitate the ecological interpretations; and keeping constant the biological 
parameters can be thought as to consider the effects of spatially-dependent 
perturbations on an averaged model of growth (TjQ), where the parameters are 
averaged in space and time. Advanced mathematical theory of averaging for 
partial differential equations can then be used to support this strategy and 
to make robust conclusions derived by such approach, at least in the case of 
small amplitude of oscillating parameters r(t,x), D(t,x) and K(t, x) (see e.g. 
Hale and Verduyn Lunel . 199Cll ; Chekroun and Rogues . 20061 ) . 



In the forthcoming computations, we assume that harvesting starts on a pre- 
viously not harvested population, which has reached its stable positive steady 
state. In other words, the environment is assumed to have reached its carrying 
capacity at t — 0, the time at which harvesting is started, with w(0,x) = K, 
where u is the solution of ([2]). This assumption is a natural one for studying 
the destabilizing effects due to removal terms. Our model (T5]) being introduced 
and discussed, the analysis of environmental effects brought by the removal 
term Y is addressed according two harvesting strategies that we present now. 



2.2 Harvesting strategies 



The first type corresponds to a "quasi" -constant-yield harvesting, with a re- 
moval term 

Y(x,u) = 6-x(x)Ps(u). (3) 
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This model is close to but different from the so-called threshold harvesting 
mod el found in ordinary d i fferen tial equations (ODEs) models of harvesting 
(see iRoques and Chekrounl . 120071 ) . 



In ©, x( x ) is a function taking the value 1 if x belongs to a harvested region, 
and the value if x belongs to a protected region. It constitutes therefore what 
we call the harvesting field, with 5 a positive constant which corresponds to 
the harvesting intensity in this field. The parameter 5 can also be interpreted 
as a quota. The last term p £ (u) is a density-dependent threshold function: 

p E (s) =0 if s < 0, p e (s) = s/e if < s < e, and p e (s) = 1 if s > e, 

where e is a small threshold below which harvesting is progressively with- 
drawn. With such a harvesting function, at each location x, the yield is con- 
stant in time whenever u(t, x) > e. Note that from a mathematical point 
of view, the function p F ensu res the nonnegativity of the solutions of ([2]) 
( Rogues and Chekroun . 20071 ). Considering constant-yield harvesting func- 
tions without this threshold value would be unrealistic since it would lead 
to harvest on zero-populations. 

The second type of strategy corres ponds to a more standard proportional 
harvesting situation (INeubertl . 120031 ). in a spatial context, where 



F(x,n) =£- X (xK (4) 

The function % is defined as above, and the term E ■ x( x ) can now be inter- 
preted as a harvesting effort at the location x. The instantaneous yield at a 
point x is then proportional to this effort and to the local population density. 



2.3 The model of fragmentation 



The re exist several ways o f obt aining hypot hetical landscape distributions, see 
e.g. iGardner et al.l (119871 ) and iKeittl (120001 ) for neutral landscape models, and 



Mandelbrotl (119821 ) for measures of fragmentation based on fractal dimension . 



The model retained here is the one developed by IRoques and Stoical (120071 ). 
and inspired from statistical physics. This is a neutral landscape model in the 
sense that it is a stochastic model of landscape pattern, and the value - pro- 
tected or harvested in the present case - assigned to a position in the pattern is 
a random variable. As key property, this stochastic model provides a numerical 
procedure for generating several samples of landscapes with breaking apart of 
the habitat while keeping constant habitat abundance, along with an exact 
control of the later and the type of breaking apart. This model offers therefore 
an appropriate framework for assessing the eff e ct of fragmentation per se on 
RD models, in the sense underlined by iFahrigl (120031 ). We make precise here 
the main parameters calibrating the model for our present purpose. 
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(a) s(x) = 94 (b) s{x) = 274 (c) s( X ) = 460 

Fig. 1. Some samples of harvesting field configurations. The black areas correspond 
to protected regions, where X = 0- 

The harvested and protected regions are entirely determined by the harvesting 
field x( x ) present in (J2J via the removal term Y . In order to build gradually 
fragmented configurations of these regions, we have discretized the domain Q 
into 50 x 50 subcells Cj, with in some cells x( x ) = f° r the unharvested 
cells, and y(x) = 1 in the other cells. Based on the stochastic model of 



Rogues and Stoical (120071 ). we have built 6000 samples of such functions x( x )j 
with different degrees of fragmentation. In all these samples, the protected 
region occupies 10 % of the domain Q. The fragmentation of the protected 
region is defined as follows. The lattice made of the cells Ci is equipped with 
a 4-neighborhood system V{Cj). We set s( X ) =number of pairs of neighbors 
(Ci,Cj) such that X takes the value on Cj and Cj. This number s(x) is 
directly linked to fragmentation: the protected region is all the more aggre- 
gated as s(x) is high, and all the more fragmented as s(x) is small (Fig. [1]). 
Therefore, s(x) can be seen as an "aggregation index" of the protected re- 
gion. On our samples, the aggregation index s varies from 94 to 460. For each 
aggregation index incremented as follows Sk '■— 94 + 6 x [k — 1), we picked 
up arbitrarily a configuration Xk with s(xk) — s k- This lead to 62 harvesting 
field distributions, with gradually aggregated configurations of the protected 
region. 



2.4 Methods 



The response to the spatial perturbation terms Y(x, u), distributed according 
to the aggregation index of their underlying harvesting fields, as described 
above, is analyzed in terms of total population size P(t), and annual yield R(t) 
in the region Q. More precisely, we evaluate the time-dependent quantities 

P(t) := / w(t,x)dx, 
corresponding to the total population at time t, and 

R(t):= C /r(x,n(r,x))dxdr, 
Jt-i Jn 

corresponding to the annual yield during the year that precedes t. 
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Note 1: In the case of quasi-constant-yield harvesting, if u(t, x) is greater or 
equal than e for all r living within the temporal window (t — l,t) and x in 
the domain Q, we simply obtain 



R(t) = 5 ■ [area of the harvested region]. 



In the case of proportional harvesting, we have: 



R(t) = -E-[mean population in the harvested region during the year (t — l,t)]. 



We estimate the intertwined effects of fragmentation per se and harvesting 
intensity by plotting P(t) and R(t) against 5 for the quasi-constant yield 
strategy, and against E for the proportional harvesting strategy. For every 
configuration, the result is represented by a curve which color is attributed 
in function of the aggregation index; red corresponds to the more aggregated 
configurations, and blue to the more fragmented configurations. For each fixed 
harvesting intensity, we computed the gap between the maximum and mini- 
mum population sizes obtained over the 62 harvesting field distributions, and 
we expressed it in terms of relative loss obtained in the worst configuration 
compared to the best one, through the formula: 100 x (highest population — 
lowest population) /(highest population). Similarly, for the annual yield, 100 x 
(highest annual yield — lowest annual yield) / (highest annual yield) was com- 
puted for each fixed harvesting intensity, where the maxima and minima are 
taken over the 62 harvesting field distributions. 



For the numerical setup we consider f2 to be a square domain of 300 km x 300 km. 
We set r = 1 year -1 and K = 10 3 individuals/km 2 . The diffusion coefficient D 
varies between 10 km 2 /year (low mobility) and 100 km 2 /year (high mobility); 
see the book of Shigesada and Kawasaki (Il997l ) for some observed values of r 
and D, for several animal species. The threshold e is set to 10 individuals/km 2 . 
Our results, except in Section 13"^ are presented at a fixed time t — 5; this time 
has been chosen to fit usual times for observations and responses in anthropic 
harvesting activities. 



The numerical integrations of the RD models were performed using a second 
order finite elements method where the solutions u(t,x) of the model (TS]) have 
been computed with the initial condition and harvesting strategies discussed 
above. The quantities of interest P(t) and Q(t) are then easily computable. 
The numerical results with ecological interpretations are discussed in the fol- 
lowing section. 
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Fig. 2. (a): Total population in Q after 5 years, in function of the quota 5. (b): Total 
yield during year 4, R(5), in function of 5. Each curve is associated with a different 
configuration of the protected region. Blue curves correspond to more fragmented 
configurations, and red curves correspond to more aggregated configurations. The 
black dotted lines indicate, for each S, the relative losses obtained in the worst 
configurations compared to the best ones. These computations were carried out for 
D = 50. 

3 Results 



3.1 Quasi- constant-yield harvesting strategy 



For every configuration, the higher the quota 5, the smaller the population 
size -P(5) (Fig. [2} a). On the other hand, a maximum yield R{5) is reached 
for an intermediate value of the quota, while small values and large values of 
S both lead to small yields (Fig. [2, b). 

As it could be expected, for small values of S, there is no dependence of the 
yield with respect to the habitat configuration. Indeed, in such cases, the pop- 
ulation density should be everywhere above e, and the harvesting function is 
therefore constant in time, equal to S in the harvested regions; see Note 1. 
However, this linear dependence of -R(5) pursues for higher values of 5 in frag- 
mented configurations. This indicates that the population density never falls 
below e on these configurations. Thus, the more fragmented configurations 
first lead to higher yields for intermediate quotas, with a yield loss in ag- 
gregated configurations which attains 88%. Conversely, the more fragmented 
configurations lead to lower yields for higher values of 5. Overall, the maximum 
yield is attained on the more fragmented configuration. Similarly, population 
sizes are higher for fragmented configurations and low quotas, with losses in 
population size up to 54% on aggregated configurations. They then become 
higher for aggregated configurations and higher quotas. It is noteworthy that, 
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whenever the harvesting term is really constant, i.e., in the region where R(5) 
is linear, fragmented configurations lead to higher population sizes. 

Such a reversion of the influence of fragmentation on both P and R, for in- 
creasing quotas, is not intuitively obvious. Yet, we can give a reasonable expla- 
nation for it. In fragmented configurations, the mean distance to a protected 
region — where population density is higher — is reduced compared to more 
aggregated configurations. Therefore, when 5 is not too large, at each location 
in Q, the population can be efficiently sustained by the protected regions, and 
the density never falls below the threshold e, leading to higher yields, com- 
pared to more aggregated configurations. For larger values of 5, the harvested 
regions become very hostile, and because of dispersion, populations tend to 
extinct, even in the protected regions. In the case of protected regions with 
small perimeters, corresponding to aggregated configurations, dispersion of the 
individuals into the harvested regions is reduced. With such configurations, 
populations can therefore sustain higher quotas without risking extinction, 
leading to higher values of P and R. Note that, whenever u(t, x) is less than 
or equal to e everywhere in the harvested region, quasi-constant-yield harvest- 
ing becomes equivalent to proportional harvesting, with effort E = 5/e. 

Comparable qualitative results were found for diffusion coefficients D ranging 
from 10 to 100 and are thus not shown. 



3.2 Proportional harvesting strategy 



When the harvesting function Y is of proportional type (j3J), our model re- 
— DV 2 u = ru(l — E ■ x(x)/r — u/ K). Population persistence 



duces to % 



for this model has been thoroughly investigated (Cantre ll and Cosnerl . 12003 



Berestycki et al. . 2005a ; Rogues and Hamel . 2007 ; Rogues and Stoica . 20071 ). 
both analytically and numerically. 

The specific effects of fragm entation of the prot e cted r egions can be deduced 
from the numerical study of iRogues and Stoical (120071 ). It shows that, on ag- 
gregated configurations, higher efforts E can be sustained without risking 
extinction. The optimal shapes of the protect ed regions, in terms of max imum 
sustainable effort, has even been obtained in (IRogues and Hamell . 120071 ). Sur- 
prisingly, these shape depend on the area of the protected region. Indeed, small 
areas have been proved to lead to disc-shaped optimal shapes, while high areas 
lead to stripe-shaped optimal shapes. 

Following our approach, we still focus on the guantities P and R, on which 
the effects of fragmentation have not yet been investigated. For every config- 
uration, the population size P(5) decreases as the effort E increases (Fig. [3j 
a). Contrarily to the guasi-constant case, with such a proportional harvesting 
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Fig. 3. (a): Total population in after 5 years, in function of the harvesting effort 
E. (b): Total yield during year 4, -R(5), in function of E. To each curve corresponds 
a different configuration of the protected regions; the blue curves correspond to 
more fragmented configurations, and the red curves correspond to more aggregated 
configurations. The black dotted lines indicate, for each value of E, the relative 
losses obtained in the worst configurations compared to the best ones. For these 
computations, we fixed D = 50. 

strategy, aggregated configurations always lead to larger populations, what- 
ever the effort. The effect of fragmentation/aggregation of the protected region 
on the population sizes really becomes noticeable when P(5) falls below one 
fourth of the environment carrying capacity. In such a case, the involved mech- 
anisms are those of persistence, which require in particular the state 0, where 
no individuals are present, to be very repulsing (o r equiv alently "unstable" , 
see lShigesada and Kawasaki! . ll997l ; lBerestycki et al.l . l2005al ). It is therefore not 



surprising to obtain effects of fragmentation comparable to those described in 
the existing literature. 

On the other hand, the yield P(5) again reaches a maximum (cf. Fig. [3j b). As 
in the quasi-constant-yield case, the maximum yield is attained for the most 
fragmented configuration. Fragmented configurations have a longer perime- 
ter, and thus provide higher transfer rates into the harvested regions. Indeed, 
obtaining higher yields with smaller total populations, for small values of E, 
implies that the size of the population situated outside the protected regions 
is higher for fragmented configurations; see Note 1. In the quasi-constant case, 
such a larger "unprotected" population would not have implied higher catches 
whenever u greater than e; this explains the qualitative difference between the 
two harvesting strategies in terms of the effects of fragmentation on P(5). 

As above discussed, higher values of the effort give a significant advantage to 
aggregated configurations, in terms of total population sizes. This translates 
into higher yields on aggregated configurations. Still, these qualitative results 
do not depend on the values of diffusion coefficients, in the selected range 
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R(5) x 10 ? R(5) 



(a) (b) 

Fig. 4. (a): Total population P(5) in 0, after 5 years, in terms of the yield dur- 
ing year 4, R(5), for the quasi-constant-yield strategy, (b): Total population P(5) 
vs yield R(5), for the proportional harvesting strategy. Blue curves correspond to 
more fragmented configurations, and red curves correspond to more aggregated con- 
figurations. For these computations, we fixed D = 50. 

[10, 100], and are not shown. 

3.3 Total population vs. annual yield 

In Fig. 01 are depicted the population sizes P(5), in terms of the yields i?(5), 
and of the level of fragmentation of the protected region for the quasi-constant- 
yield and proportional harvesting strategies (Figs.Hl (a) and (b), respectively). 
The flip shapes of these diagrams teach us that for a given yield two branch 
of disjoint intervals of population size are admissible, sufficiently far to the 
left from the "bending point" . The upper branch corresponds to low harvest- 
ing intensities (below that leading to maximum yield), and the lower branch 
corresponds to higher intensities. Remarkably, for each given population size, 
higher yields are obtained on more fragmented configurations; the lower the 
population size, the higher this effect. 

This D-shape of the P — R diagram supports the idea that it will be difficult 
to predict the quantitative effect of fragmentation in practice, at a fixed yield, 
without knowledge of the total population size; a situation which typically 
arises in ecological application where such a knowledge is difficult to achieve. 

3.4 Results for times t G (0, 20) 

Fig. [5] depicts how the total population P(t) depends on the harvesting inten- 
sity and on the aggregation index in function of the time t. 
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Fragmented configurations Aggregated configurations 




(a) (b) 

Fig. 5. (a): Total population P(t) in function of time t and in terms of the quota 5, 
for the quasi-constant-yield strategy, (b): Total population P(t) in function of time 
t and in terms of the effort E, for the proportional harvesting strategy. For these 
computations, we again fixed D = 50. 



For the quasi-constant-yield harvesting strategy (Fig. a), we observe that 
the higher the quota, the sooner the inversion of the effects of fragmentation. 
Thus, the threshold quotas, above which aggregated configurations lead to 
higher populations, decrease as t increases. 



For the proportional harvesting strategy (Fig.0 b), at each time t, aggregated 
configurations are still associated with larger populations. Moreover, the effect 
of fragmentation tends to increase with time. 



In both cases, the amplitude of the effect of fragmentation depends on the 
harvesting intensity and on the time point of the analysis. However, those 
results, and in particular the quasi-constant-yield case, suggest previous sec- 
tions results might still be qualitatively true at times other that t — 5. This is 
confirmed by numerical computations (not presented here) which indeed lead 
to comparable results, but with a shift in the harvesting intensity. 



Analytical studies show that the solution u( t, x) of ([2]) converges to some equi- 
librium state for both quasi-constant- yield (IRoques and Chekrounl . 120071 ) and 
proportional ( Berestycki et al. . 2005al ) harvesting strategies; however, conver- 
gence rates are not known. Figs. [5] (a) and (b) provide information about these 
convergence rates: the higher the harvesting intensities, the sooner equilibria 
are reached. Thus, at t = 5, population sizes may be almost at equilibrium for 
high harvesting intensities or still farer to reach it for lower intensities. 
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4 Concluding remarks 



In the existing single-species modeling literature, fragmentation of the favor- 
able region is usually found to be detrimental for population survival. On 
the other hand, recent empirical studies show that fragmentation may have 
positive or negative effects, and that positive effects occur more often. 



In her review paper on the topic iFahrigi (120031 ) proposed two reasons for 



negative effects of fragmentation of these favorable regions. Firstly, favor- 
able patches are too small to sustain a local population, and secondly, the 
total perimeter of the favorable region is large, leading to increased transfer 
rates into the unfavorable regions. She also proposed several reasons for posi- 
tive effects of fragmentation. Among these reasons, sm aller distance bet ween 



patches, higher immigration rates into the patches (see iGrez et all 120041 ) and 
positive edge effects for some species are evoked. 

In our work, based on single-species reaction-diffusion models with harvesting 
terms, we have captured both positive and negative effects of fragmentation 
of the favorable region, interpreted here as a protected region. 

Firstly, for large harvesting terms, we found that aggregated configurations of 
the protected region lead to higher population sizes, and higher yields than 
fragmented configurations. In that sense, our results are not contradictory 
with previous modeling results: under hostile conditions, the extinction risks 
are higher, and the chances of persistence are increased on aggregated config- 
urations. The reasons evoked by Fahrig for negative effects of fragmentation 
then become paramount in our model at low population sizes. 

On the other hand, fragmented protected regions lead to higher population 
sizes when a constant number of individuals are removed per unit of time 
in the unprotected region (constant-yield-harvesting). Furthermore, for fixed 
population sizes, the two harvesting strategies we studied in this paper lead 
to higher yields for fragmented configurations; cf. Fig. HI Such higher yields 
may stimulate harvesters to slow down, and may therefore be beneficial for 
populations. 



Sumailal (119981 ). through a compartment modeling approach, showed that 
transfer rates between favorable and unfavorable regions were of critical im- 
portance for understanding the role of reserves, and the author emphasized the 
necessit y of a more precise mod eling approach of the transfer rate function. 



Latter, iTischendorf et al.l (120051 ) found, via a simulation approach, negative 
effects of fragmentation when the probability of individuals to go from the 
favorable to the unfavorable region was high, and positive effects in the oppo- 
site case with a high unfavorable to favorable boundary- crossing probability. In 
our models, the diffusion coefficient is spatially-constant, and the boundary- 
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crossing probabilities are therefore equal in inward and outward directions. 
However, net transfer rates from favorable to unfavorable regions depend on 
the geometry of the protected regions, and on the relative population densities 
inside and outside the favorable regions. Indeed, in dividuals mov es are driven 
by random diffusion, and Green's formula (see e.g. Evanj . 1998 ) implies that 



the instantaneous population flux from the protected region to the harvested 
one is: 



^-Protected ^harvested region = ~D j y 



du , 

—ds, 

on 



where T denotes the boundary of the protected region, and |^ is the outward 
gradient in population density observed on this boundary. Our work shows 
that both the geometry of the protected regions, and the relative population 
densities inside and outside the favorable regions interact to give negative 
effects of fragmentation for high harvesting terms, i.e. when the contrast be- 
tween protected and harvested areas is strong, and positive effects for less 
contrasted environments. 

In multi-species models, fragmentation p e r se c an alter interactions among 



species. As reviewed by iRyall and Fahrigi (120061 ) . predator-prey models pre 



diet varying effects of fragmentation on equilibrium densities of predator and 
prey populations, depending on the specific assumptions of these models. In 
this paper we have demonstrated that single-species reaction-diffusion models 
with removal terms can support relative effects of fragmentation per se. The 
methods and results could serve as a first step to bridge the gap between em- 
pirical work and modeling, for react ion- diffusion as well as other single-species 
models. 
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